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Histograms of two images of the fundus. 



(57) Abstract: The present invention comprises a 
robust technique , to automatically grade retinal images 
through detection of lesions that occur early, in the 
course of diabetic retinopathy: dot hemorrhages or 
microaneurysms, blot and striate hemorrhages, lipid 
exudates and nerve-fiber layer infarcts. In addition, 
the present invention includes methods to extract the 
optic nerve in the appropriately identified fields, and 
to track and identify the vessels (measuring vessel 
diameters, tortuosity and branching angles). The present 
invention preferably identifies 3 levels: no retinopathy, 
microaneurysms alone, and lesions additional to 
microaneurysms, the latter two levels being the earliest 
detectable forms of the disease. The method of the 
present invention is particularly suited to overcoming 
the difficulties in grading retinopathy, which stems 
from image-to-image variation, low contrast of some the 
lesions from the background and non-uniform lighting 
and flare within the same image resulting in variation in 
the different quadrants of the same image. The system 
has the capability of utilizing human expert marking 
of images and grading of retinoathy to improve image 
quality assessment, lesion detection and retinopathy 
grading. 
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ance Notes on Codes and Abbreviations" appearing at the begin- 
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System and Method for Screening Patients for Diabetic Retinopathy 

Summary of the Invention 

5 

One of the complications of diabetes is the gradual degradation of vision known as 
diabetic retinopathy. It is known that over time, a majority of diabetics will lose some 
vision to this condition, and the present state of the medical art is. to treat the retinal 
lesions that mark the disease with laser light. Diagnosing diabetic retinopathy,, 

10 however, requires a trained specialist to view the retina (or a photograph of the retina) 
at periodic examinations, and to recognize small lesions and changes therein. 
Because there are many more diabetics than specialists, and because they must be 
examined at regular intervals, difficulties and deficiencies have arisen in health-care 
systems that attempt to screen and treat diabetic retinopathy. The present invention 

15 provides an efficient computer-implemented screening analysis of retinal photographs 
to identify the stages of diabetic retinopathy. 

Background of the Invention 

20 Diagnosis of diabetic retinopathy is commonly performed by skilled medical 
personnel (usually ophthalmologists) in direct examination of the fundus, or by 
evaluation of sets of fundus photographs taken using special-purpose cameras. Such 
examinations, whether done directly or by review of photographs is a time-consuming 
and inexact process, and experts in diagnosis often disagree in their results, 

25 particularly in the very earliest stages of the disease. In addition such diagnosis is 
expensive, and is often performed at intervals that are longer than desired due to the 
cost and lack of available skilled manpower. 



30 



Brief Description of the Invention 

The present invention comprises a robust technique to automatically grade retinal 
images through detection of lesions that occur early in the course of diabetic 
retinopathy: dot hemorrhages or microaneurysms, blot and striate hemorrhages, and 
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lipid exudates. In addition, the present invention includes methods to detect nerve 
fiber layer infarcts, extract the optic nerve in the appropriately identified fields, and to 
track and identify the vessels (measuring vessel lumen diameters, tortuosity, and 
branching angles). The present invention preferably identifies 3 levels: no 
5 retinopathy, microaneurysms alone, and lesions additional to microaneurysms, the 
latter two levels being the earliest detectable forms of the disease. The method, 
however, may be expanded to 7 levels through the detection of the lesions that occur 
in advanced stages of the disease and may be utilized, through the evaluation of 
changes of the normal vasculature (lumen diameter, tortuosity, and branching angle), 
10 to evaluate the risk of developing retinopathy. 

The method of the present invention is particularly suited to overcoming the 
difficulties in grading retinopathy, which stems from image-to-image variation, low 
contrast of some the lesions from the background and non-uniform lighting and flare 
15 within the same image resulting in variation in the different quadrants of the same 
image. 

The expert decision system implemented in the present invention bases the 
retinopathy grade determination upon the results of the lower-level detectors (lesion 
20 detection) for each eye photographed. The system may be tuned for individual 
detectors and for mydriatic images and non- mydriatic images using separate 
^ parameters based on the camera type, image type, characteristics of patients, and 

characteristics of fundus. 

25 A data archive is the core of the data management system that acts as the central 

repository for all patient data: images, demographics as well as reports. This archive 
will be accessible for storage and retrieval via the Internet and will be of a scalable 
design to accommodate growth over time. The benefits of a centralized data 
management architecture include (1) The ophthalmologist or retina specialist will 

30 have access to current as well as prior studies for comparison regardless of where the 
historical data was acquired. (2) A central data repository will allow for an objective 
and quantitative means to evaluate the progression of the disease in individuals or 
populations over time in terms of the changes of normal vessels, the change of 



WO 03/020112 



PCT/US02/27586 



existing diabetic lesions ,and the occurrence of new ones. This database will provide 
the foundation to develop regression and serial studies to develop risk prediction 
algorithms in the future. (3) Algorithms that will scan the archive will produce 
quantitative measures of vascular tortuosity, branching angle, and caliber variation, 

5 which have been identified as markers of vascular disease and which can be tracked 
over time. Again, this will also enhance risk prediction, predominantly in the early 
stages, prior to or subsequent tothe development of retinopathy. It should be noted 
that these parameters are not able to be assessed by human grading. (4) Data mining 
of the massive warehouse of data can allow screening proficiency and patient 

10 compliance to be examined, as well as providing valuable insight into the trends in the 
various populations and a comparison to observe and follow over time the 
effectiveness of treatment among patient populations (e.g. similar patients in each 
clinic). 

15 Brief Description of the Figures 

Figure 1 depicts histograms of two images of the fundus. 

Figure 2 depicts an initial fundus image, and the result of retina extraction processing. 
Figure 3 depicts the result of processing the retina image by coalescent filter [Bhasin 
20 and Meystel 94] to obtain the vessels. 

Figure 4 depicts graphically the problem with locating the retina in a photograph. 
Figure 5 depicts insetting. 

Figure 6 (a) depicts a subsampled original retinal image. 

Figure 6 (b) depicts the image of Figure 6(a) after subtraction of vessels from the 
25 image. 

Figure 7 depicts results of finding the retina using geometry to improve the pickup. 
Figure 8 depicts an overall flow diagram for the method of the present invention. 
Figure 9 depicts a flow diagram of the image processing method of the present 
invention. 

30 Figure 10 depicts a flow diagram of the photographer assistance method of the present 
invention. 

Figure 1 1 depicts a flow diagram of the serial (over time) study method of the present 
invention. 
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Detailed Description of the Invention 

According to the present invention, images are acquired locally with immediate 
feedback to photographer to assist in optimizing each image quality using a small 
looping process to guide the photographer, analyzed locally and stored locally on CD- 
ROM and magnetic storage. Screening is reported immediately from the system with 
a printed report. Those eyes that "fall out" of screening, either because a significant 
number of photographs cannot be analyzed (at present too poor despite 3 attempts to 
photograph) for one eye or because they are above the accepted threshold in 
number/type of lesions/level of retinopathy, will be transmitted to a remote specialist 
ophthalmologist for reviewing. The specialist will then transmit back a review with 
grading and recommendations: either the patient is to be advised to make an 
appointment for further examination, or the photographs are sufficient but all the eye 
needs is repeat screening at a recommended interval. This report is maintained in the 
database for data mining with regard to the success and thresholds for screening and 
referral. 

In practicing the present invention, images may be obtained in one of the following 
20 ways: 

(a) Monochromatic 30° images digitized from 35 mm Ektachrome® slides at at 
least 1 024 X 1 024 X 8-bit depth using a green Wratten filter 

(b) Photograph eyes with a digital on-line system using a fundus camera with a 
. digital camera back allowing color photography at at least 2 meg x 2 megx 32 bit 

25 depth resolution 

(c) Photograph eyes with a digital on-line system using a fundus camera with one 
or more interposed filters including 535nm notch filter, 605 nm notch filter (or 
others), a digital camera back allowing monochromatic 1 1024 X 1024 X 8-bit depth 

30 The present invention preferably uses a protocol of 5 photographs of 35° fields per 
eye. 




-4- 



WO 03/020112 



PCT/US02/27586 



The present invention also provides for cataloging patients (including naming images 
as well as providing other demographic and clinical data) to enable batch processing 
and historical tracking of disease as well as identification of risk indicators. It has 
been determined that cataloging patients is important: keeping digital photos on file is 
efficient for doctor's office - for example, images may be maintained in an image- 
management database which also allows the recording of ancillary patient data and 
physician notes for each set of photgraphs. 

The image-management system of the present invention allows medical professionals 
to track improvement and determine more accurately what factors in what patient 
classes dictate degradation rate, which dictates patient recall and reexamination 
frequency. 

The system of the present invention uses a file-naming convention of the images that 
allows the system to automatically locate images of the same patient and the 
corresponding field. The naming scheme for all images is: PatientID- Camera 
type.Eye-Field-Field of view.?Image type.Processing 

• Eye is either L or R 

• Field is a number from 1 thru 5 

• PatientID is system generated alphanumeric identifier 

• Camera type is a number 1 for non-mydriatic camera, 2 for mydriatic camera 

• Field of view is the angular degree of retina photographed by camera with each 
image (i.e. 30° or 35° or 45 °) 

• Image type is a number for the type of image: 1 for digital direct color image, 2 
for monochromatic derived from digitization from Ektachrome slides, 3 for 
monochromatic followed by the notch filter peak: (e.g. 535 for a 535nm peak 
filter). 

Processing is a code that indicates the type of the image: RAW (for original 

photograph), GRA (graded image indicating lesions present). 

The method of the present invention includes a software system for assisting the 

photographer in capturing adequate retinal images. In use, the photographer 

photographs the eye and indicates the field being photographed. When these images 

are captured, each photograph is marked with an identifier: the current system utilized 
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but not to be restricted i.e. OD field#l centered on fovea with the disc on the right or 
left of the field depending up on the eye; OD field #2: disc in lower left (superonasal) 
OD field #3: disc in upper left (inferonasal); OD field #4 disc in upper right corner, 
fovea in upper center (inferotemporal; and OD field #5 disc in lower right corner, 
5 fovea in lower center (superotemporal); the OS fields are just the opposite, right from 
left. These 5 photographs will give some overlap between the 45 degree images with 
several photographs taken of the fovea and span all of the necessary areas of the retina 
that require photogrpahy and examination. It is known that the present invention may 
be improved if the photographer marks the fovea in field #1 . Each photograph is 
10 then: 

• Examined for these necessary elements, disc and vessels and their position in the 
photograph etc. 

• Contrast: The background is compared with the vessels (loss of contrast is 
produced by moving the camera too far in or too far out- on the Canon CR6, for 

15 example, this is helped by haying the photographer align some bright alignment 

dots. . . . ■ 

• Focus: de-focused image is detected by examining edge definition of central 
vessels. 

• Alignment (Up/down/right/left within the pupil): Misalignment with regard to the 
20 pupil produces an edge flare within the field of the photograph (loss of contrast):- 

this will require contrast evaluation at the edge of the field as well as in the center 
to detect the flare (to be compared with a database of some example images, taken 
but yet to be determined) and tell the photographer to realign if the edge flare is 
detected. 
25 ! 

If the photograph is determined to be inadequate for any of the parameters examined, 
a set of instructions for each element is immediately presented to the photographer 
with instructions about which field and eye to re-photograph, and how to improve the 
image. This leads, not only to improved quality of images acquired for an individual 
30 but also leads in general to more rapid improvements in retinal photography by the 
inexperienced photogapher. 
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Retinopathy Grading 

Retinopathy grading in the present invention is based on image processing and 
computer vision software for lesion detection and quantification from digitized 
monochromatic images taken either on line with 1280X1024 pixel density digital 
monochromatic camera, or utilizing the green channel from a trichromatic color 
camera with similar one channel pixel density, or digitized at similar resolution from 
slide film images at 30° or 50° taken from a non-mydriatic camera. Reference is made 
to the CD-ROM Appendix which contains the computer source code of the preferred 
embodiment of the present invention. 

Quantification: (for each field) number /field, total area/ field area evaluated, 
histogram of number vs size, histogram of number vs. density in each field. 



Screening results in 3 levels: 1) no lesions (hemorrhages, lipid exudates, or nerve 
15 fiber layer infarcts), 2) microaneurysms only (dot hemorrhages), 3) dot hemorrhages 
with other hemorrhages, exudates or infarcts. 

The separate detectors are for the following lesions and anatomical features in the 
retina: 

20 A. Hemorrhages: dot, blot, striate 

B. Lipid Exudates 

C. Nerve fiber layer infarcts 

D. Optic Nerve Head 

E. Vessels (Arteries, Veins): Detection of 1°,2°,3 0 branching order 

25 

Other lesions of more advanced retinopathy grades for which detectors are being 
pursued: 

A. Intra-retinal microvascular anomalies (IRMA) 

B. Venous loops 

30 C. Epi-papillary neovascularization 

D. Epi-retinal neovascularization 

E. Sub-hyaloid or vitreal hemorrhage 

F. Epi-retinal fibrosis 
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G. Retinal detachment 

Approaches to determine grade include: 

1 . Geometric techniques 

2. Image Processing to obtain Position, Size of Lesions 

3. Matched Filter - Find something that is oval - use Biology code 

4. Morphological Filter - Grow-Shrink method 

5. Remove vessels - 

-Find all junctions 

- Find vessels emanating from them and follow them 

Region of Interest Extraction 
For images digitized from film, the portion containing the retina is extracted. For 
both digital and digitized images, the bad contrast portions caused due to flare, etc. 
are removed. An improvement to the present invention removes certain portions of 
photographs which result from image flare and contain no useful data. 

The histograms of two different fundus images shown in Figure 1 demonstrate that 

there is a considerable separation in gray values between the background and the 

* - 

retina. Unfortunately this separation occurs in a broad range of values.. There are 
other factors such as flare at the edges of the retina and artifacts in the imaging that 
severely limit the use of a simple technique. 

These vessel boundaries will be used to determine the primary, secondary and tertiary 
vessels. An adaptive threshold is used to locate the vessel. Its settings depends on 
image contrast. This requires locating the retina within the image and measuring its 
contrast. To discriminate the retina from the background the inventors have used the 
ISODATA clustering technique to first define a threshold value and then binarized at 
that value. The ISODATA algorithm is an iterative method based upon the gray level 

histogram of the image. The histogram is split up into two parts, the foreground 

» 

pixels and the background pixels, assuming an initial threshold value. Then the 
average value of the foreground and of the background pixels is calculated and a new 
threshold value is taken midway between these two average values. This process is 
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repeated, based upon the new threshold estimate until the threshold value does not 
change any more. 

It often happens that this does not result in definition of the full retina because of the 
5 presence of a deep shadow over part of the image (Figure 4). One way to overcome 
this obstacle is to use a restricted convex hull. First the largest object is found and 
then all the combinations of two Contour points with a Euclidean distance less than or 
equal to some preset distance (d), are connected by a straight line. If a background 
pixel is found on such a line, it is added to the original object. This operation also 
10 closes all holes which are less than d pixels wide. 

An alternative approach is not to use the convex hull, which creates the need to close 
off any internal holes within the segmented retina. To do so, generate a binary image 
using the ISODATA technique and remove all but the largest object. Then invert the 
15 binary image and remove all objects that are not contiguous with the boundary of the 
image, that is that are not part of the background. Then re-invert the binary image. 

Once the retina is located, inset the retina mask edge to define an inner boundary (Fig 
5). An adequate inset distance for the lk images is 100 pixels. The reason for 
20 performing this operation is to avoid including the retinal edge in measuring the retina 
contrast Often the edge is subject to deep shadow or bright flare.. 

The histogram of the retina mask following insetting is used to define the region of 
the image for which the gray value histogram is determined. From the histogram the 

25 standard deviation is computed and that is used as a measure of the image contrast. 
The rank order to be utilized in the obtaining of the vessels needed to be derived. 
Different ranks were used and the results are shown for rank = 0.25, 0.5 (median 
filter) and 0.75 in Figure 6. Here the original image was a subsampled version of the 
scanned image to reduce the processing time. The subsampling was by a factor of 4. 

30 The subsampling was done by a straightforward linear method and other techniques 
like bilinear sampling might improve the results and will be studied later. To see the 
effectiveness of the technique, the vessels found were subtracted from the original 
image (subsampled scanned image by a factor of 4). The original image and the 
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results of the subtraction are presented in Figure 7. We see that there are a lot of 
false alarms which can be isolated by form factor and total area of object. 

Discrimination of the vascular bed from the background is accomplished on the basis 
5 of density, size and shape The algorithm utilizes an adaptive threshold based on a 
rank order filter. The routine assigns a value of 0 (background) or 1 (object) to a 
given pixel depending on its value relative to that of its neighborhood: All the pixels 
of the neighborhood are rank ordered on the basis of gray value, a constant offset is 
added to the value of the central pixel and the result compared to the value of the fa 
10 element of the neighborhood. If it is lower, that is the central pixel is substantially 
darker than its surround, then the central pixel is assigned a value of 1, else it is given 
a value of 0. Three interrelated parameters important to the operation of this adaptive 
threshold are the offset, kernel size and count threshold (i ). The kernel size and offset 
are constants both for a given image and across images. We have tested both circular 
15 and rectangular kernels. As you would predict, the circular ones result in less artifact. 
We settled on a circular filter with a diameter of 17 pixels. The offset used was 15. 
The count threshold is kept constant for a given image but is varied from image to 
image depending on the image contrast. 

20 To set the count threshold the standard deviation of the retinal image is used. To 

assess how the count threshold varies with the image contrast, the standard deviation 
of the histogram of each of twenty-five images was computed and each image was 
convolved with the adaptive threshold with the count threshold varied until the best 
segmentation was achieved as judged visually. The relation between the two has been 

25 modeled by a linear fit (y=234.34-1.607x; r = 0.88; Fig. 3). 

Vessels identified and removed from image to leave only dot hemorrhages as dark 
structures. The present invention can remove objects like the optic nerve head which 
are very bright to get a clearer representation of the gray value distribution 
30 (histogram) in the rest of the image. Much of the processing is dependent on grey 
values and these techniques along with specialized algorithms like the k-means, 
modified gradient, signature analysis and morphology allow us to perform adaptive 
segmentation that leads to a more reliable extraction of the lesions. 
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The present invention applies a smoothing filter - special configuration that removes 
noise but does not destroy the structure of the hemorrhages or cause artifacts in the 
image. 

The lesions detected by the various detectors are fed into an expert system of known 
construction. The expert system is constructed out of rules built from consultations 
with expert ophthalmologists. It embodies the reasoning used by them to make 
interpretations. Among the features queried of the experts and rule-coded are: 



1 . Quick scan of the image that detects of vessel fragments. 

2. Morphological processing to connect co-linear fragments and branches to the 
major vessel (dilate and then constrict the regions) 

3. Skeletonization algorithm developed for finding the median line 

15 4. Subtract the vessels found from the original image to see efficacy of processing 



The present invention uses a path-searching technique. To find the path, we start at 
the position in the gray value image corresponding to a given endpoint and searche for 
the direction in which the mean gray value is maximum. The initial search direction is 

20 determined by examining the skeleton and the connectivity of the endpoint and setting 
the initial direction (d) to the opposite direction. Given a current pixel (x, y) and a 
current direction (d), there are three candidate points for the next pixel, viz. the 
neighbor points in the directions (d-1), (d) and (d+1) module 8. To decide which to 
take, three straight lines are considered, one in the direction (d), another in between 

25 the directions (d) and (d+1) and a third in between the directions (d) and (d-1). Along 
each of these lines the mean value is calculated, over a length of <val> pixels. The 
next pixel is the candidate point that corresponds with the direction with the 
maximum average. The process is repeated at the new pixel, using the found direction 
as a new value for (d). The process stops, firstly, if the distance of the candidate pixel 

30 to the image boundaries is less than a preset value. Further, if the maximum of the 
mean values in the candidate directions is less than the specified value <minval>. 
Next, if a pixel is obtained that is set in the output bitplane. For this the system uses 
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the background outside the eye. Finally the process is stopped if the number of pixels 
found becomes equal to a preset value. 

Alternatively examination of the gradient normal to the generated curve at a regular 
5 spacing may be used to reject the curve. 

Gradient 

CStatistics stats(m_doc, original); 
CBasicFilter bf(m_doc, m_src); bf.Sobel(); 
10 CStatistics stats2(m_doc, m_src); 

BYTEthr = (BYTE)(stats2.m_mu + kl* stats2.m_sd); - 
bf.Threshold( thr ); 

Copy either m_src->m_data[i] or 0 to modest depending on 
contrast measure of background (to eliminate bright flare) 
15 if( ref->m_data[i] > (stats.m_mu + k2* stats.m_sd) ) 
m_dest->m_data[i] = 0; 

else 

m_dest->m_data[i] = m_src->m_data[i]; 
20 Methodology for Detectors 

■ 

Separate the veins from arteries, measure vessel branching patterns, vessel tortuosity, 
vessel calliber variation. Locate the disc and fovea so that lesions can be localized in 
position relative to the fovea. 



25 



30 



General philosophy: Encode tests as functions that return three possible scores: 

1 (SUCCESS), 0 (Indeterminate or DONT_KNOW) or -1 (FAIL) 
This setup will be used later to multiply output of functions with weights and combine 
to produce a confidence score. 

As discussed above, following the adaptive threshold and noise removal, vessel 
segmentation is achieved on basis of object size and form factor. The objects rejected 
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are in fact classified as hemorrhages. Size is used to distinguish dot from blot 
hemorrhages (DH & BH, respectively) with later being the larger of the two. 

Adaptive Segmentation 
5 Run K-means on large image segments - 4x4 or 8x8 segments in whole image - to 
' produce three regions: background, dark regions like vessels and dot hemes, and light 
regions like lipid exudates. This will give a rough segmentation of the image that 
takes the local shading and variation into account. The mean grey values of the 
regions will be used later as pivotal values for finer segmentation. 

0 

Technique for hemes 
If any 2 of following tests are failed, i.e. if (test 1 0 + test2() + ... == -2) 

0. Rank order the grey values in a NxN region and then divide them into M-tiles (if 
M = 5, it would be quintiles). The pixels in the lowest rank order are colored 

L5 RED. ' / 

• For a dot heme, the RED part should be somewhat in the center as the ideal 
profile of the heme is a crater - low grey values in the center surrounded by 
concentric rings of larger value. 

• A lipid exudate is like a mound which has the same characteristics as heme for 
20 an inverted image. 

1. Centre 4x4 of full stamp does not have any RED 

2. Let XMEST, YMIN, XMAX, YMAX be corners of full stamp and (ymin, xminl, 
xmaxl) be first ival of RED, and (ymax, xmin2, xmax2) be last ival - these form 
actual comers of RED. Check if any edges of RED are at comers of full as in the 

25 case for triangles formed in peripheral segments of image and in flare regions As 
you have a regionarray, XMIN = lesionpost->supp.xl, ... so test is 

(lesionpost->supp.xl - xminl > -2 and lesionpost->supp.yl - ymin > -2) OR 
(lesionpost->supp.xl - xmin2 > -2 and lesionpost->supp.y2 - ymax < 2) OR 
(lesionpost->supp.x2 - xmaxl < 2 and lesionpost->supp.yl - ymin > -2) OR 

30 (lesionpost->supp.x2 - xmax2 < 2 and lesionpost->supp.y2 - ymax < 2) 

3. To check if it is a vessel, see if it runs for full length of stamp 

lesionpost->supp.yl - ymin > -2 and lesionpost->supp.y2 - ymax < 2 

4. If total size of stamp is > 60x60 return FAIL . 
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5. Form factor of RED should be between 0.9 and 1.1 which eliminates trapeziods 
for vessels and triangles and allows circular objects 

Technique for Cotton Wool Infarcts and Lipid Exudates 

The CWS and EX are differentiated from the background of the fundus images by 
first inverting the image gray value so that the exudates are now darker than the 
background and then applying the same adaptive filter used to segment the 
hemorrhages. After removal of noise with a median filter the remaining objects are 
separated into CWS and EX. This is done on the basis of size and gradient. Objects 
smaller than a fixed constant and with a sharper gradient are classified as EX. 
Remaining objects are classified as CWS. Obviously, the gradient is assessed on the 
original gray value image. 

Refinements 

Developed an approach to detect objects which are touching the vessels. It employs a 
distance transformation which replaces each pixel of an object with an estimate of its 
shortest distance to the background (the distance to the nearest background pixel). 
To speed production of the map, the inventors have used pseudo-Eu'clidean distances: 
distance to the nearest neighbor is 5, to the diagonal it is 7, and knight move is 1 1 . 
Now thresholding at some value, for instance, 50, will capture pixels which are 10 
pixels or further from the nearest edge of the object. In this way it is possible to 
detect an unusually thick vessel. 

A refinement of this approach looks not at the absolute thickness but rather the change 
in thickness. Here we are looking for a sudden change from a thin to a thick vessel 
and then back down to a thin vessel. To do this the skeleton of the binary vessel 
image is found using a distance transform. It differs from the medial-axis transform 
in that here the skeleton is defined as a set of connected, one pixel thick arcs, lying 
midway between the object boundaries and being a topological retraction with the 
same connectedness as the original object. Thus unlike the medial axis connectivity is 
preserved. It differs from the Hilditch skeleton in that the way it is applied, the 
distance skeleton is computed using pseudo-Euclidean distances rather than city-block 
distances (1 nearest neighbor, 2 diagonal, 3 knight move). Now, the skeleton is used 
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as a mask with the distance map to obtain distances from the skeleton to the edge of 
the object. 

Level 1 

5 Filter and Run CLesion (Dot,BIot etc) with 

IndicationFilter<DarkonLight/LightonDark> 

(Converts grayscale to binary with small blobs in image space). 

Filter to smooth noise and generates vessel bed. The use of filtering is configurable 
10 and is determined through the use of experiments done on images from a specific 
camera that measure the relative merit. 

Runs a filter with 3x3 cross at center and 8-neighbor 'arms* ending in 2 pixel long 
'hands' at end of each arm. The center grey values (weighted sum) is compared with 
15 sum of 'hands'. If darker than all, hit is marked at center 3x3. Hands are moved out by 
StepSize from center, starting at minSize and ending at. Multiple passes were added 
to pick up more than one lesion size; maxSize. Is crude matched filter - many;false 
alarms. Call with minSize & maxSize ODD, StepSize EVEN. Looks like: 



20 



* * * 

* * * 



* * * 

25 * * * 



Last step is to PruneBySize... 

Converts CRa to object list and applies simple geometrical rules to prune objects in 
list. List is held and managed in CIArray - see definitions there, min and max 
30 threshold are used to specify color-range to consider for object e.g. red and green is 
from 1 to 2. Allows pruning by: # pixels constituting object, and min & max extent of 
object. Permits specification of AND/OR for x & y extent-limit rule and whether or 
not to FILL objects that are passed. 
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Level 2 

Processes in object space. Initiates call to <CObjetsManipulation::PruneBySize>. 
PruneBySize is a common function for all lesions, has no separate body 

BOOL CLesion : :GeometricFiIter(BOOL bUseDlg) 

Converts CRa to a cohesive group of pixels or a list object and apply simple 
geometrical rules to prune the objects in the image by: # pixels constituting object, 
and min & max extent of object. 

Compute first and second Moments(area, centroid, major and minor axis of best-fit 
ellipse) 

Testl : area > minPix 

Test2: Compare both xExtent & yExtent to minimum and max thresholds or check if 
only either of them meet the minimum constraints. Permits specification of AND/OR 
on the constraint i.e. either dimension or both dimension to be constrained for x & y 
extent-limit rule. 

Level 3 

Refines image-space and object-space evidence with help of expert rules. 

CLesion::SignatureFilter(BOOL bUseDlg) 

Modified 02/28/97 to allow processing based on cf5 OR raw signature. Major 
overhaul of CDialogManager also done, and ini file also. 3 new stringsections per 
lesion ADDED to ini and 1 DELETED i.e. L3[lesion]Dlg replaced by 
L3[lesion]CommonDlg, L3 [lesion] Cf5Dlg, and L3[lesion]RawDlg. First Dig 
determines whether user wants to use cf5 or raw (BOOL, m_paraml0) and other 
parameters common to both methods. Second dialog is selected on the basis of value 
of m_paraml0 and queries for those parameters which are specific to the particular 
method selected earlier. 

CLesion: :SignatureFilterBasedOnRaw 
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Function SPLIT 2/27/97 to allow different processing of two types - original (using 
cf5) or new (using raw) for analysis Function modified 1/17/97 to create diagnostic 
file with fixed name in current results directory. This file contains yellow filled, 
passed objects for final analysis. File is overwritten on each pass through. Valid 
5 during batch processing of a single image. 

PROCESS 

int reason; 
BOOL bPassed; 

10 CScreenerApp* pApp = (CScreenerApp*)AfxGetAppO; 

CRegionarray* refCRa = new CRegionarray(m_doc, m_doc->m_pFiles.raw) 
CRegionarray* vslCRa = new CRegionarray(m_doc, m_doc->m_pFiles.vsl) 
CRegionarray* diagCRa = new CRegionarray(m_doc, refCRa->m_dX, refCRa- 
>m_dY) 

15 ' 

CIList ci(m_doc) . . , 

CIArray ca(m_doc, m_doc->m_pImage, minThr, maxThr)Get objects formerly in 

"bigHst" 

20 while( ca.GetNext( &ci ) ) 

. { 

BOOL bNotOn Vessel = NewNotOnVessel(&ci, vslCRa, &bPassed, &reason); 1 Pixel 
overlap 

if( bNotOnVessel || (m_lesionType=lBLOT) ) 
25 { 

CExpertSystem ce(m_doc, m_lesionType, bNotOnVessel, &ci, refCRa, userbins, 
borderWidth, bUsePercentile, nObjectsMax, filffactorMin, formFactorMin, 
fillFormMin, pcntAreaMin, trl 1, tr21, tr31, trl2, tr22, tr32, bUseLocalContrast, 
minGValLT, maxGValGT, minDRange, maxDRange, minBinWidth 1 , 
30 minBinWidth2, minBinWidth3, minBinWidth4); 



if( (bPassed = ce.Verify(&reason)) = TRUE ) 
{ 
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InsertLesionInOutputImage(&ci, diagCRa, ce.m_analysis->m_j>tamp, 
bStampBorderOnly, bStampOneColor); 

ce.UpdateDocStatisticsO;For each incremental PASSED lesion 
5 record info for entry in database 

} 

WriteOutputlmage(diagCRa); 



10 



30 



Default constructor creates an array of SIGNATURESTATISTICS that is currently 
HARDCODED to size = 5; Of these index 0 is a placeholder for statistics of the 
15 stamp's CRegionarray parameters and 1-4 holds statistics for the R/G/Y/B stamps 

(islands in regionarray). In verification, it will be necessary to guarantee that statistics 
of each color have meaning. NULL structs in the array will be checked there. 

CExpertSystem: :CExpertSystem(CScreenerDoc* pDoc, int lesionType, BOOL 
20 notonvsl, 

CIList* ci, CRegionarray* refCRa, int userBins, int borderWidth, BOOL 
bUsePercentile, int nObjectsMax, double fillFactorMin, double formFactorMin, 
double fillFormMin, double pcntAreaMin, int trl 1, int tr21, int tr31, int trl2, 
int tr22, int tr32, BOOL bUseLocalContrast, int minGValLT, int maxGValGT, 
25 int minDRange, int maxDRange, int minBinWidth 1 , int minBinWidth2, int 
minBinWidth3, int minBinWidth4) 

1 . Instantiate CStamp Aanalysis (which looks at <RED if dot/blot, YELLOW 
lipid/cwool (OIS), BLUE lipid/cwool (Kaiser)> part of signature) 



2. Initialize stats 

typedef struct _SIGNATURESTATS 
{ 
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shortnObjects;# of objects (of a particular color) in a given stamp (CIArray) 
doubleaveSize;average size 
doublesizeSD;modal size 

doubledispersionjcentroid of areas normalized by center x,y 
doublelargestArea;pixels in largest object 
CRectlargestRect;rect of largest object 
} SIGNATURESTATS ; 

SIGNATURESTATS m_ss ? SetSize(5);Array of SignatureStatistics structures 

Initialize:m_ss[0].nObjects, .largestArea, .aveSize, .sizeSD, .dispersion, 
.largestRect.left, .largestRectright, .largestRecttop, .largestRect.bottom 
(This corresponds with stamp Rect, the other numbers are not used) 



15 GetSignatureStatistics . 

Function added 030497 to assist in designing an expert system for signature analysis. 
SIGNATURESTATS definition is in CScreenerApp. Parameters to be added till 
accuracy of the damned system is adequate - or unimprovable. May split into multiple 
hierarchical functions if too much redundant processing becomes necessary with a 

20 single function. Will keep a single SIGNATURESTATS structure in any case. 

m _ss[ 1 ] = ((CIArray*)(m_analysis->m 

((m_analysis->m_stamp->m__x 1 +m_analysis->m_stamp->m_x2)/2), 
(m_analysis->m_stamp->m_yl+m_analysis->m_stamp->m_y2)/2); DO RED 



25 



30 



3. Do Signature Analysis Including 

Investigates 8-neighborhood of core-color object to count and record the transitions to 
next-color bands. In other words looks for deviation from ideal case which is 

44444444 
43333334 
43222234 
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43211234 
43211234 
43222234 
43333334 
5 44444444 

This cannot be done with list objects because "next-to" and "nearly-surrounds" are 
hard to compute. This function provides a working metric on a scale of zero to 8 for 
transitions from each band into the next, e.g. (7,6,3) which allows use in rules. 
10 Function requires start point within core-color region of regionarray. Yet to determine 
effect of concave objects. 

FAIL unconditionally if local contrast not high enough and sanjay wants to use this 

rule Added 070998 
15 if( m_bUseLocalContrast ) rule to apply test or skip it 

if( m_stamp->m_bwl <= m_minBinWidthl ) 
{m_bSuccess = FALSE; m_reason = 22; return;} 

20 if( (m_stamp->m_gvMax - m_stamp->m_gvMin + 1) < m_minDRange ) 
{m_bSuccess = FALSE; m_reason = 23; return;} 

if( (m_stamp->m_gvMax - m_stamp->m_gvMin + 1) > m_maxDRange ) 
(m_bSuccess = FALSE; mjreason = 23; return;} 
25 } 

FAIL unconditionally if largest object bleeds 

m_bOnEdge = CheckEdge(m_ss[ 1 ] .largestRect, m_ss[0] .largestRect); 
if( m bOnEdge ) {m_bSuccess = FALSE; m_reason = 6; return;} 

30 

FAIL unconditionally largest object too small 

m_prcntArea = m_ss[l].largestArea/ max(1.0, m_ss[0].largestArea); 
if(m_prcntArea < m_pcntAreaMin ) { m_bSuccess = FALSE; m_reason = 5; return;} 
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FAIL unconditionally if too much clutter 

if(mj5s[l].nObjects > mjtiObjectsMax ) { m_bSuccess = FALSE; m_reason = 1; 
return;} 

FAIL unconditionally if shape skewed 



m_dYl=m_ss[l].largestRect.bottom-m_ss[l].largestRect.top; 
mJbrmFactor = min(m_dXl, m_dYl) / max(l .0, max(m_dXl, m_dYl)); 
mJillFactor = m_ss[l]. largest Area / max(1.0, m^dXl * m_dYl); 
if(m_formFactor < m_formFactorMin ) { m_bSuccess = FALSE; m_reason = 2; 
return;} 

if(m_fillFactor < m_fillFactorMin ) { m_bSuccess = FALSE; m_reason = 3; return;} 

SUCCEED if InvestigateNeighborhood comes up with a "good" {a,b,c} gradient 
int xStart = ( m_ss[l].largestRect.right +m_ss[l].largestRect.left )/2; 
intyStart==(m_ss[l].largestRect.bottom+m_ss[l].largestRect.top )/2; 
InvestigateNeighborhood(m_stamp, xStart, yStart, FALSE, &m_trCount 1 , 
&m_trCount2, 
&m_trCount3); 

if(m_trCountl >=m_trCountlMinl ) 
if( mJrCount2 >= m_trCount2Minl ) 
if( m_trCount3 >= m_trCount3Minl ) 
{ m_bSuccess = TRUE; mjreason =10; return; } 

i£( rnJrCountl >= m_trCountlMin2 ) 
if( m_trCount2 >= m_trCount2Min2 ) 
\{{ m_trCount3 >= rnJrCount3Min2 ) 
{ m_bSuccess = TRUE; m_reason =12; return; } 

4. After each passed lesion 



m_dXl= m_ss[l].largestRect.right -m_j>s[l].largestRect.left; 



void CExpertSystem: :UpdateDocStatisticsQ 
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switch( m_lesionType ) 
case 1DOT: 

(m_doc->m_fs. nDots)++ ; 

(m_doc->m_fs.aveDotDefl2) += (double)m_trCountl; 
(m_doc->m_fs.aveDotDef23) += (double)m_trCount2; 
(m_doc->m_fs.aveDotDef34) += (double)m_trCount3; 
break; 
Etc. 

5. After each processed image (all lesions) write to mdb for postprocessing (Grading) 

BOOL CDiagnosis::UpdateFieldMDB() 
CFldResults pFR; 

pFR-AddNewO; 
pFR.m_FieldCode = m_fCode; 
pFR.m_sdGrayVal - m_doc->m_fs.sd; 
pFR.m_meanGrayVal = m_doc->m_fs.mu; 
pFR.m_minGrayVal - (BYTE)(m_doc->m_fs.minG); 
pFR.m_maxGrayVal = (BYTE)(m^doc->m_fs.maxG); 
pFR.mnDots = m_doc->m_fs.nDots; 
pFR.m_nBlots = m_doc->m_fs.nBlots; 

pFR.m_nLipids = m_doc->m_fs.nLipids; - 
pFR.m_nCwools = m_doc->m_fs.nCWools; 
pFR.m_aveDotDefinition = 

100*(int)(m_doc->m_fs.aveDotDefl2/max(1.0,m_doc->m_fs.nDots)+0.5)+ 

10*(int)(m_doc->m_fs.aveDotDef23/max(1.0,m_doc->m_fs.nDots)+0.5)+ 

(int)(m_doc->m_fs.aveDotDef34/max(1.0,m_doc->m_fs.nDots)+0.5); 

pFR.m_aveBlotDefinition = 1 00*(int)(m_doc->m_fs.aveBlotDefl 2/max(l .0,m_doc- 
>m_fs.nBlots)+0.5)+ 

10*(int)(m_doc->m_fs.aveB16tDei23/max(1..0,ni_doc->m_fs.nBlots)+0.5)+ 
(int)(m_doc->m_fs.aveBlotDef34/max(l .0 > m_doc->m_fs.nBlots)+0.5); 
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pFR.m_aveLipidDefmition = 100*(int)(m_doc- 

>m_fs.aveLipidDefl2/max(1.0,m_doc->m_fs.nLipids)+0.5)+ 

1 0*(int)(m_doc->m_fs.aveLipidDef23/max(l .0,m_doc->m_fs.nLipids)+0.5)+ 

(int)(m_doc->mJs.aveLipidDeO 

pFR.m_aveCwoolDefinition = 100*(int)(m_doc- 

>m_fs.aveCWoolDefl2/max(1.0,m_doc->m_fs.nCWools)+0.5)+ 

1 0*(int)(m_doc->m_fs.aveCWoolDef23/max(l .0,m_doc->m_fs.nCWools)+0.5)+ 

(int)(m_doc->m_fs.aveC W6olDef34/max( 1 .0,m_doc->m_fs.nCWools)+0.5); 

pFR.UpdateO; 

6. After mdb is updated for patient batch, all field results from the database (7/eye) 
are unified to generate a grade between 1 & 3. 



Quantify severity and use of Neural Net to generalize 

15 Neural Net improves segmentation based on features extracted from the image 

* Histogram of number of lesions and size of lesions 

* Types of lesions 

The best kind of matched filter is a classical neural network. 
20 [a] For each lesion, pick an input layer sized to match the lesion size 

[b] Train on the areas in the target 

[c] Apply to the candidate regions identified during the false positive biased test. 

[d] There is no conceivable reason to identify more than one kind of lesion in a single 
pass, or to view the image as a single entity. 

25 [e] The problem of classifying existence and number of lesions is fundamentally 
different from your dissertation because we don't need to recognize a "holistic" 
relationship or an "entity", just small objects with "clearly" defined properties in the 
bitmap. (The big objects are not relevant at this stage of the work, we can probably 
gloss over the features which require cross referencing amongst themselves). 



30 



ENTITY (dynamic programming based segmentation) modified to accept a penalty 
function for the area of image to be parcelled); intializing penalty function (show with 
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colors in the annular region before starting PARCEL); defining prohibition zones (set 
penalty =1) 

Extraction of lesions given image and indication of field 
Intelligent agent 

Algorithm that learns params from images marked by a retina specialist. Associates 
params with image type. 

This will help in having in one JAZZ, the following for 10 patients: the original, 
expert (dot/blot, lipid, striate, cotton wool spot), 4 comparisons and intermediate 
results 

Implementation issues 
Several polymorphic data structures representing regions and lines such that different 
representations enable us to write efficient algorithms for different kinds of 
processing. For example, a run-length encoded representation of a binary image 
allows fast determination of statistical properties like area, whereas an array 
representation allows one to access random cells faster for morphological processing. 
Image Encoding Method for efficient lookup for recognition and interpretation 

Physical implementation of the system is envisioned in the form of: 

I) Database of elementary images in SQL. 

II) Access through SQL to the Expert Validated Knowledge Base. 

III) Database of scripts for finding each region in LISP: Rules for all Types of 
Knowledge. 

IV) Inference Engine in LISP: OTTO, OPS5, etc. 

V) Low level routines in C: Filters; ENTITY (dynamic programming based 
segmentation) modified to accept a penalty function for the area of image to be 
parcelled); intializing penalty function (show with colors in the annular region before 
starting PARCEL); defining prohibition zones (set penalty = 1) 
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Comparison of sequential images 

Comparison of sequential images to improve the risk prediction for pathologic 
changes occurring over time is provided by overlay and comparison of feature 
differences and by comparison in the database of the number and location of lesions 



for instance toward the fovea). 

While the present invention has been described above in terms of specific 
embodiments, it is to be understood that the invention is not limited to the disclosed 
10 embodiments. On the contrary, the present invention is intended to cover various 
modifications and equivalent structures included within the spirit and scope of the 
appended claims. 



5 



(e.g. to detect new lesions or their migration toward anatomical features of interes : 
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We claim as our invention: 

1 . A method for diagnosing diabetic retinopathy by computer implemented 
screening of retinal photographs comprising the steps of: 

(A) receiving one or more images of a human retina and storing the retinal 
images in the memory of a digital computer; 

(B) processing the stored original images to identify the presence of features 
indicative of diabetic retinopathy, the features of selected from the group of dot 
hemorrhages, microaneurysms, blot hemorrhages, striate hemorrhages, and lipid 
exudates; and 

(C) reporting the absence or presence and nature of the features indicative of 
diabetic retinopathy. 

2. The method of claim 1 further comprising selecting those stored images for 
processing that are deemed acceptable in terms of at least one image quality criterion 
selected from the group of discernible necessary elements, contrast, focus, alignment, 
and completeness. 

3. The method of claim 1 wherein the received retinal images are captured 
from a human retina by a retinal camera. 

4. The method of claim 3 further comprising: 

(D) assessing the stored retinal images in terms of at least one image quality 
criterion selected from the group of discernible necessary elements, contrast, focus, 
alignment, and completeness; and 

(E) if the assessment is unacceptable in any of the selected image quality 
criteria of step (D), then prompting an operator to capture further retinal image. 
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5. The method of claim 1 where in the retinal photographs are taken over one 
or more predetermined intervals and compared across one or more of the intervals to 
assess the development of the features indicative of diabetic retinopathy. 

5 6. The method of claim 1 where in the original photographs are stored in the 

memory of the digital computer using an indexing convention that denotes the identity 
of the patient, the eye and field imaged, and the processing applied. 

7. The method of claim 1 wherein one or more images are selected based on 
10 predetermined criteria and transmitted for examination by a human expert. 



8. The method of claim 7 wherein the human expert is provided with both the 
transmitted images and the reported features indicative of retinopathy from step (C) of 
claim 1. 

15 . . ' ^ • 

9. The method of claim 7 wherein the result of examination by the human 
expert is transmitted using the network and is stored in relation to the transmitted 
images. 

20 10. A system for diagnosing diabetic retinopathy comprising a computer 

including a processor and memory wherein the processor is programmed to 

(A) receive one or more images of a human retina and storing the retinal 
images in the computer memory; 

(B) process the stored original images to identify the presence of features 
25 indicative of diabetic retinopathy, the features of selected from the group of dot 

hemorrhages, microaneurysms, blot hemorrhages, striate hemorrhages, and lipid 
exudates; and 

(C) report the absence or presence and nature of the features indicative of 
diabetic retinopathy. 
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1 1 . The system of claim 10 further comprising a retinal camera for capturing 
images from a human retina which is operatively linked to the computer. 

5 12. The method of claim 10 further comprising at least one additional 

computer operatively linked to the computer using a network. 



13. The method of claim 12 wherein one or more stored retinal images are 
selected based on predetermined criteria for transmission over the network for 
10 examination by a human expert. 



14, The method of claim 13 wherein the human expert is provided with both 
the transmitted images and the reported features indicative of retinopathy from step 
(C) of claim 10. 



15 



15. The method of claim 13 wherein the result of examination by the human 
expert is transmitted using the network and is stored in relation to the transmitted 
images. 



-28- 



WO 03/020112 ^ PCT/US02/27586 

1/11 




BEST AVAILABLE COPY 



PCTAJS02/27586 

2/11 




Figure 2 . (a) Initial image; (b) result of retina extraction 
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Figure 4. Problem with locating retina. 



BEST AVAILABLE COP 



* 



WO 03/020112 



PCT/US02/27586 



5/11 




background 



inset contour 



Figure 5. Insetting. 
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